% -----------
% defining the variables and parameters
% -----------
close all;
var y pai r rr z g epsr ym paim rm;

varexo etaz etag etar;

parameters
beta paistar rstar tau kappa phii phipai phiy old_ylead old_ylag old_rpai
new_ylead new_ylag new_rpai
rhoz rhog rhor
lambda elas mu teta corinov gammap gammac;

 
% -----------
% setting the parameter values
% -----------



% Log data density is -364.324934.
%  
% parameters
%            prior mean post. mean conf. interval prior pstdev
% 
%       lambda   0.300   0.1794  0.1027  0.2471 beta  0.1000
%       gammac   0.500   0.4997  0.3920  0.6205 beta  0.2000
%         elas   3.000   2.7236  1.9729  3.4768 gamm  0.5000
%       phipai   1.500   1.8578  1.5800  2.1311 gamm  0.2000
%         phii   0.500   0.6452  0.5654  0.7308 beta  0.2000
%         phiy   0.250   0.1168  0.0248  0.2092 gamm  0.1500
%      paistar   4.000   3.9006  2.5337  5.3430 gamm  2.0000
%        rstar   2.000   2.8814  2.3330  3.4660 gamm  1.0000
%         rhog   0.700   0.8631  0.7938  0.9248 beta  0.1000
%         rhoz   0.700   0.7491  0.6639  0.8347 beta  0.1000
%  
% standard deviation of shocks
%            prior mean post. mean conf. interval prior pstdev
% 
%         etag   0.380   0.2021  0.1505  0.2486 invg  0.2000
%         etaz   1.000   0.8780  0.7391  1.0246 invg  0.5200
%         etar   0.310   0.2779  0.2358  0.3229 invg  0.1600
%  
% correlation of shocks
%            prior mean post. mean conf. interval prior pstdev
% 
%    etag,etaz   0.000   0.5781  0.4240  0.7468 norm  0.4000

%Steady State Values
lambda=0.1794;
elas=2.7236;
teta=0.75;
rstar=2.88;
beta=(1+rstar/100)^(-1/4);
paistar=3.9006;
mu=0.2;
gammap=0.5;
gammac=0.4997;
%Policy Parameters
phii=0.6452;
phipai=1.8578;
phiy=0.11680;


%AR(1)
rhog=0.8631;
rhoz=0.7491;

%Standard errors

sdz=0.8780;
sdg=0.2021;
sdr=0.2779;
corinov=0.5781;


% -----------

model(linear);

y=((1-(lambda*elas/((1-lambda)*(1+mu)))*(1+gammac*mu/(1+elas*(1-gammac))))
/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*y(+1)
+((gammac*(1-lambda*elas/((1-lambda)*(1+elas*(1-gammac)))))
/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*y(-1)
-((1-gammac)/(1+gammac-(lambda*elas/((1-lambda)*(1+mu)))*(1+(gammac*(1+2*mu))/(1+elas*(1-gammac)))))*(r-pai(+1))+g;

pai=(gammap/(1+(1+rstar/100)^(-1/4)*gammap))*pai(-1)
+((1+rstar/100)^(-1/4)/(1+(1+rstar/100)^(-1/4)*gammap))*pai(+1)+
((1-teta)*(1-(1+rstar/100)^(-1/4)*teta)/teta)*((1/(1-gammac))+(elas/(1+mu)))*(y-z)
-((1-teta)*(1-(1+rstar/100)^(-1/4)*teta)/teta)*(gammac/(1-gammac))*(y(-1));


r=phii*r(-1)+(1-phii)*(phipai*pai(+1)+phiy*y)+epsr;

z=rhoz*z(-1)+etaz;

g=rhog*g(-1)+etag;

epsr=etar;

ym=y;

rr=r-pai;

paim=paistar+4*pai;

rm=rstar+paistar+4*r;

end;

  


steady;

check;
shocks;
var etag = sdg^2;
var etaz = sdz^2;
var etar = sdr^2;
var etag,etaz = (corinov*sdz*sdg);
end;

periods=1000000;
stoch_simul(order=1,irf=0,simul);

